High velocity spikes in Gowdy spacetimes 
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We study the behavior of spiky features in Gowdy spacetimes. Spikes with velocity initially high 
are, generally, driven to low velocity. Let n be any integer greater than or equal to 1. If the initial 
velocity of an upward pointing spike is between 4n — 3 and 4n — 1 the spike persists with final 
velocity between 1 and 2, while if the initial velocity is between 4n — 1 and 4n + 1, the spiky feature 
eventually disappears. For downward pointing spikes the analogous rule is that spikes with initial 
velocity between 4n — 4 and 4rz — 2 persist with final velocity between and 1, while spikes with 
- initial velocity between An — 2 and 4n eventually disappear. 
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There have been several investigations of the approach to the singularity in inhomogeneous cosmologies, predomi- 
nantly in the case that a two-dimensional symmetry group acts spatially [jj, |2j, L&Jji Lj- Lj- III- La- L3- USi llil LL2l • The most 
extensively studied class of such spacetimes is the class of Gowdy spacetimes 13] on T 3 x R. Numerical studies of 
these spacetimes show that the approach to the singularity is asymptotically velocity term dominated (AVTD) except 
at an isolated set of points. 

As shown in [(| this behavior can be understood by considering certain terms in the equations as "potentials" 
that affect the dynamics at each spatial point. These potentials drive the dynamics into the AVTD regime, except 
at those isolated points where one of the potentials vanishes. The behavior at each of these isolated points differs 
from that of its neighbors, leading to the creation of features called spikes that become narrow as the singularity is 
approached. The analysis of uses an explicit approximation to find a closed form expression for the behavior of 
q I the spikes. However, this approximation only remains valid throughout the approach to the singularity for a certain 
O" 1 class of spikes: those for which the "velocity" v of the spike is "low" (1 < v < 2 for "upward pointing spikes" and 
< v < 1 for "downward pointing spikes"). What then happens if a spike forms and it initially has "high" velocity? 
bi). Certain terms in the evolution equations which are decaying for low velocity spikes are instead growing in magnitude 
for a high velocity spike and the sign of the net effect is such that the velocity should decrease 0, . An explicit 
approximation that is valid as the velocity at a high velocity spike decreases has not been found, leaving open the 
■ question, "What is the long term behavior of initially high velocity spikes in Gowdy spacetimes?" 

To answer this question, we perform numerical simulations of Gowdy spacetimes with these high velocity spikes. 
- - ■ Since the spikes become ever narrower as the singularity is approached, until just this narrowness causes the velocity 
to decrease if it is originally high, extremely high resolution is needed to follow the detailed behavior of the spike. We 
use a numerical method that was used in | 14| to study critical collapse. Here the equations are put in characteristic 
form and the outermost grid point is chosen to be the ingoing light ray that hits the singularity at the center of the 
spike. In this way, the grid gets smaller as the feature that it needs to resolve shrinks. 

Section 2 presents the equations and the characteristic numerical method. Results are given in section 3 and 
conclusions in section 4. 

II. EQUATIONS AND NUMERICAL METHODS 

The Gowdy metric on T 3 x R takes the form 

ds 2 = e^H-^i-dt 2 + dx 2 ) + t[e p (dy + Qdz) 2 + e - p dz 2 } (1) 
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where P, Q and A are functions of t and x. The vacuum Einstein field equations split into "evolution" equations for 
P and Q 

Ptt + r l P, t - P xx + e 2P (Q 2 x - Q 2 t ) = (2) 

Q.u + t' 1 Q. t - Q,xx + 2{P,tQ,t - P X Q, X ) = (3) 

and "constraint" equations for A 

A, f = t[P 2 + P 2 X + e 2P (Q 2 t + Q 2 X )\ (4) 

\ x = 2t(P x P t + e 2P Q, x Q,t) (5) 

(here jfl = d/da). The constraint equations determine A once P and Q are known. The integrability conditions for 
the constraint equations are satisfied as a consequence of the evolution equations. Since the evolution equations do 
not depend on A there is essentially a complete decoupling of constraints from evolution equations. Therefore, for the 
purposes of this paper we will treat only equations (|2I3|) . The only restriction that the constraints place on initial 
data for equations ill'IMI) is the following: since A at x — is the same as A at x = 2ir, it follows that the integral from 
to 2-7T of the right hand side of equation JSJ must vanish. We require that this restriction is satisfied by the initial 
data for equations ()2l3ll and then these equations insure 1 that the restriction is also satisfied at subsequent times. 

The singularity is at t = 0. It is often helpful to introduce the coordinate r = — lni. Thus the singularity is 
approached as r — > oo. In terms of this coordinate, the evolution equations (I2I3II become 

P.tt ~ e 2P Q 2 T e- 2 ^P xx + e 2 ^Q 2 x = (6) 
Q,tt + 2P T Q, T - e- 2r (Q. xx + 2P X Q. X ) = (7) 



The velocity is defined as v = \jP 2 r + e 2P Q 2 T . In the limit of large r, if we neglect the terms in equations H6I7|) 

proportional to e~ 2r then we find that v — > Voo and P — > v^t where Woo depends only on x. However, if Voo > 1 then 
neglect of the term e 2 ^ p ^ T ^Q 2 x is not justified. What happens is that if P T is greater than 1, then this term drives it 
to values less than 1, except at the isolated points where Q. x — 0. At those points upward pointing spikes form. The 
previous studies [a,Sll3 suggest that if 1 < P T < 2 at the spike, the spike will persist, but that if P T > 2 at the 
spike, P TT will become negative and its magnitude will start growing exponentially in r as the spike narrows. But 
the previous studies leave open the question, "What happens next?" We want to know the outcome. 

For our purposes, we will need the field equations in characteristic form. To that end we introduce null coordinates 
u = t + x and w = t — x and characteristic variables A, B, C and D given by 

A = (u + w)P, w (8) 

B = {u + w)P, u (9) 

C=(u + w)e p Q, w (10) 

D = (it + w)e p Q, u (11) 

The evolution equations (12l.'it then become 

A, u = (u + w y 1 [CD + (A-B)/2] (12) 
B iW = (u + w)~ 1 [CD + (B-A)/2] (13) 

C, u = (u + wy^-AD + iC -D)/2] (14) 

D, w = (u + wy 1 [-BC+(D-C)/2] (15) 

We impose the condition that P(t,x) and Q(t,x) are even functions of x. This insures that Q, x {t, 0) = and 
therefore that the spike will be at x = 0. With this condition, we only need to evolve on the domain x > 0. This 
corresponds to u > w. The characteristic initial value formulation for this system is as follows: on a u — constant 
surface we give as initial data A and C as functions of w. Equations (|13fl and (|15|) are then integrated to yield B 
and D on the initial data surface. The numerical method used is second order Runge-Kutta. The initial condition 
for this integration is that B — A and D — C at u — w. This is a consequence of the fact that P and Q are even 
functions of x. Now with A, B, C and D known on the initial surface, equations i|12|) and (|14f> can be regarded as 
ODEs for each value of w. At each grid point, these ODEs are integrated to yield A and C at the next value of u. 
The numerical method used is a second order predictor-corrector method. This entire process is iterated to produce 
the full evolution. 
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FIG. 1: spike in P at r = 9 

We choose the minimum value of w to be zero, corresponding to the light ray that will hit the singularity (f = 0) 
at the center of the spike (x = 0). The maximum value of w is the current value of it, corresponding to x = 0. Thus 
as the evolution proceeds, some grid points leave the physical domain and are no longer part of the evolution process. 
Throughout, we choose du — dw. This means that we lose one grid point at each time step. When half of the grid 
points are lost, we put them back in between the remaining ones and interpolate A and C to obtain their values on 
the new grid points. Thus we make the grid twice as fine. In this way, as the singularity is approached we always 
have enough grid points to resolve the ever narrowing spike. 



III. RESULTS 



While our characteristic method (c-code) allows arbitrary initial data for A and C (as long as the initial data does 
not determine a Cauchy surface 1 ), it is helpful to make comparisons with codes using the usual Cauchy methods. We 
therefore also use codes that evolve equations (12131) (t-code) and equations (|6I7[) ( T-code). These codes use standard 
centered differences for spatial derivatives and the iterated Crank-Nicholson method [l5j (ICN) for time evolution. 
We begin with Cauchy data at t = ir and evolve using the t-code to t = 7r/2. In the course of this evolution, we obtain 
A and C on the null surface given by u = 7r/2 and < w < 7r/2, which we use as initial data for the c-code. 

To obtain solutions with high velocity spikes, we note that in the polarized case (Q = 0) the velocity can take on 
any values. We therefore choose initial data with small Q. The form of the data is P = 0, P. T = scosx, Q — qcosx 
and Q :T = where s and q are constants. Figures ^ an( i [21 show the results of evolving these data. Here s = —8, q = 
1.0 x 10~ 3 and the result is produced by using the T-code and evolving to r = 9. 

Now we take the same initial data and evolve it for a longer time, first using the t-code to generate the appropriate 
initial data for the c-code and then using the c-code to evolve very close to the singularity. The result for P T vs r at 
x = is shown in figure|31 As the spike forms, P T > 2. As the spike narrows, P T decreases, as predicted. Note that 
P T at the spike starts out larger than 2, and that it is eventually driven below 2. 

One might wonder what mechanism is responsible for driving the spike velocity below 2. The answer comes from 
an examination of equation (jSJ. Rewriting this equation with P TT on one side we have 



p 



-2t 



P. a 



„2P 



c,2(P- 



(16) 



1 If a Cauchy surface is determined, then the initial data must be compatible with the periodicity of A. 
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FIG. 2: spike in P >T at r = 9 




FIG. 3: P >T at the spike 



At the spike the term e 2 ( p ~ T 'Q 2 x vanishes and the term e 2P Q 2 T is positive definite. The term e~ 2r P xx is negative, 
since the spike in P is upward pointing. What is happening is that the spike is becoming so narrow and thus P\ xx 
is becoming so large that despite the smallness of e~ 2r the quantity e~ 2T P tXX is not negligible. The term e 2P Q\ 
is also not negligible, but its magnitude is less than the magnitude of e~ 2r P iXX) so P TT is negative. In the explicit 
approximation, e 2P Q\ and e~ 2r P jXX & r e decaying in magnitude if P T < 2, but if P T > 2 and if Q. xx does not vanish 
at the spike they are growing in magnitude, and the leading part of e 2P Q 2 T is half the magnitude of the leading part 



of e 2t P iXX , so the net effect is a decrease in P TT From figure|3|we see that this decrease eventually comes 
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FIG. 4: P, r at the spike for different initial values of P :T . 

to an end while P r is still greater than 1. Now the conditions for the explicit approximation of the spike to be valid 
throughout the approach to the singularity are satisfied, indicating that the spike will persist. 

Another family of initial data which we use to generate spikes is t = to, P = p, P T = s, Q = qcosx and Q tT = 
where to, p, s and q are constants, with s > 1. The constants are chosen so that the explicit approximation to the 
spike is initially valid, until the spike narrows enough that the velocity begins to slide down. To allow matching of the 
t-code to the c-code, the value of to is chosen to be Ann/m, with n + 1 the number of grid points for the c-code and 
the m + 1st grid point identified with the 1st grid point in the t-code (corresponding to a difference in x of 2tt). In 
figures FTTH we show results from a family of initial data sets with n = 160 and m = 2050 (which results in to as 0.98), 
p = and q = 0.01. We vary s, from 2.4 to 20.7. In figure |1| we again plot P r vs t. In figure[5]we plot both P T vs r 
and the velocity, v vs t for a single initial data set, with s = 13.2. In figure|Hlwe plot v vs t. When 2 < wo < 4 there 
is just one "slide" in the velocity, v. Its final value is about 4 minus its initial value, and the final value of P r is about 
equal to 4 — s. When 4 < s < 6, there is a slide, after which P r as 4 — s (i.e., it is negative), and v w |4 — s\. After 
this slide the velocity is approximately constant, while P T "bounces" to a value equal in magnitude to 4 — s, but 
now positive (|4 — s|). This sort of bounce has been discussed in previous work and is driven by the term e 2P Q 2 T 
in the evolution equation for P. When s > 6, there is not just one slide in P T , but a number of slides. Each slide 
approximately satisfies the rule v — > |4 — v\ and P T — > 4 — P T . If P T is negative after a slide, then there is a bounce 
as just described above. If P T > 2 after this bounce, then a slide follows, and so on. In each case, finally, < P T < 2, 
P r « « and the conditions are satisfied such that the explicit approximation of should be valid throughout the 
remaining approach to the singularity. These results and all other results we obtained from various choices of initial 
data 2 agree with the picture we have just presented. 

The rule for the final (after all the slides and bounces) value of P r , based on the initial value of P r is shown in 
table U with n > an integer, and < a < 1. Both initially and finally, v « P r , but this relation does not hold 
during bounces and slides. The rule for the final value of P T based on the initial value of P r when Q x ^ 0, obtained 
from the explicit approximations Q (if P T > 1, P r — > 2 — P T and then if P T < 0, P T — * — P T and reiterate) is also 
given in the table. 3 

To obtain a picture of what happens at downward pointing spikes, one could use the same combination t-code 



2 This method allows study of a large number of different spikes with varied initial conditions into the "asymptotic regime" , since the 
code runs very quickly. 

3 The discussion in this paper assumes "genericity conditions", for example, that Q tXX ^ at upward pointing spikes. For discussion of 
persistence of high velocity spikes under special conditions, see 0. 
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FIG. 5: P :T and the velocity at the spike for P :T — 13.2 initially. 
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FIG. 6: Velocity at the spike for different initial values P T . 



and c-code with the appropriate initial data. But this is not necessary, since an upward pointing spike is mapped to 
a downward pointing spike with p^_ own = 1 — p^p by the Gowdy to Ernst transformation 0- Thus for downward 
pointing spikes an upward slide in P T occurs if P >T < —1 at the spike, and the rule for the slide is P T — > — P T — 2. 
If P T > 2 after the slide then there is a downward bounce P T — > 2 — P T . Reiterating this until — 1 < P r < 1, we 
obtain the rule shown in table ITT1 with n > and < a < 1. The Gowdy to Ernst transformation maps the condition 
Q.x = for an upward pointing spike to form to the condition e 2P Q jT = for a downward pointing spike to form. 
The c-code only gives us information concerning an extremely small neighborhood of the spike. This is precisely 
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TABLE I: Final P r from initial P jT > 0, for integer n > and < a < 1. 



Initial P )T 


Final P T if Q, x = 
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TABLE II: Final P T from initial P jT < 0, for integer n > and < a < 1. 



Initial P T 


Final P T if Q, T = 


Final P T if Q, r / 


— An — a 


-a 


cr 


—An — 1 — a 


-1+cr 


1 - a 


-An - 2 - a 


cr 


cr 


-An - 3 - a 


1 - cr 


1 - a 



its advantage. But the rules just obtained raise the question, "How does the small neighborhood of the spike match 
onto the rest of space when the rule for the final value of P r at the spike agrees with the rule for the final value of 
P T elsewhere (i.e., in half the cases)?" To see how a spike sits in the rest of the spacetime, we use Cauchy evolution 
using the symplectic PDE solver [|J with mesh refinement We use a resolution of 1024 (the 1025th grid point is 
identified with the 1st, corresponding to a difference in x of 27r). Two levels of mesh refinement are made at spikes, 
each increasing the resolution five-fold. This code is not as fast as the combination t-code and c-code, and is not as 
effective over as large of a range of initial data sets nor as far into the asymptotic regime as the combination t-code 
and c-code. But all the results obtained suggest that when the rule obtained for final value of P r at the spike matches 
the rule obtained elsewhere, the signatures of the spike are effectively washed out, for both upward and downward 
pointing spikes. We show an example, for a downward pointing spike, in figure [7] The initial data is n = 162 and 
rn = 2048 (so i ~ 0.99), s = 5.7, p = 0, q — 0.01 and now Q = qcos(x + tt/5) so that no spike is at the edge of 
the grid, for easier implementation of the mesh refinement. P T « —3.7 when the downward pointing spike begins 
to form. The spiky features are pronounced in the middle row of figure but by the time of the final row they are 
disappearing. 



IV. CONCLUSIONS 



We have seen that there is a mechanism that, in general, eventually drives upward pointing spike velocities below 2 
and downward pointing spike velocities below 1. The general asymptotic behavior of Gowdy spacetimes is then that 
given in [||. The spacetime is everywhere AVTD. At a set of isolated points there are spikes. The closer a point is to 
a spike, the longer it takes to reach the VTD regime. Upward pointing spikes in general have 1 < v < 2, downward 
pointing spikes in general have < v < 1 (—1 < P r < 0) and the shape of both kinds is described by expressions in 
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